#! /usr/bin/python
"""
evaluate in read level
"""
import argparse
import os
import random
from collections import namedtuple

import numpy as np
import gzip

# num_sites = [10000, 100000, 200000, 1000000000]
num_sites = [100000, ]

def sample_sites_ont(filenames, probcf):
    skip_cnt = 0
    cnt = {}
    group_dict = {}
    for filename in filenames:
        if filename.endswith(".gz"):
            infile = gzip.open(filename, 'rt')
        else:
            infile = open(filename)
        for line in infile:
            
            words = line.strip().split("\t")
            sampid = "\t".join([words[0], words[1], words[2], words[3], words[4], words[5]])
            
            prob0, prob1 = float(words[6]), float(words[7])
            key = tuple(words[:3]+[words[4]])  # 使用前三列+read id作为键
            if key in cnt:
                cnt[key] += 1
            else:
                cnt[key] = 1
            
            
            if key not in group_dict:
                group_dict[key] = [prob0, prob1, 1]  # [第7列总和, 第8列总和, 计数]
            else:
                #cnt[key] -= 1
                group_dict[key][0] += prob0  # 累加第7列的值
                group_dict[key][1] += prob1  # 累加第8列的值
                group_dict[key][2] += 1  # 计数加1
            
        infile.close()
    result = 0
    # 对字典中的每个键值对进行处理，计算平均值
    for key, value in group_dict.items():
        prob0 = value[0] / value[2]  # 计算第7列的平均值
        prob1 = value[1] / value[2]  # 计算第8列的平均值
        if abs(prob0 - prob1) < probcf:
            skip_cnt+=1
            continue
        
        #if prob1>prob0:
        #    label=1
        #else:
        #    label=0

        result+=1
                
    print('there are {} cpg candidates totally, {} cpgs kept, {} cpgs left'.format(len(cnt.keys()),
                                                                                   result,
                                                                                   skip_cnt))   

if __name__ == '__main__':
    parser = argparse.ArgumentParser(description='Calculate call accuracy stats of ccsmeth for cpgs')
    parser.add_argument('--input', type=str, required=True, action="append")
    parser.add_argument('--prob_cf', action="append", required=True,
                        help='this is to remove ambiguous calls. '
                             'if abs(prob1-prob0)>=prob_cf, then we use the call. e.g., proc_cf=0 '
                             'means use all calls. range [0, 1].')
    args = parser.parse_args()
    prob_cfs = args.prob_cf
    for prob_cf in prob_cfs:
        sample_sites_ont(args.input, float(prob_cf))
    